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The possibility to draw links between the isospin properties of nuclei and the structure of compact 
stars is a stimulating perspective. In order to pursue this objective on a sound basis, the correlations 
from which such links can be deduced have to be carefully checked against model dependence. Using 
a variety of nuclear effective models and a microscopic approach, we study the relation between the 
predictions of a given model and those of a Taylor density development of the corresponding equation 
of state: this establishes to what extent a limited set of phenomenological constraints can determine 
the core-crust transition properties. From a correlation analysis, we show that (a) the transition 
density p t is mainly correlated with the symmetry energy slope L, (b) the proton fraction Y Pt t 
with the symmetry energy and symmetry energy slope (J,L) defined at saturation density, or, even 
, better, with the same quantities defined at p = 0.1 fm~ 3 , and (c) the transition pressure Pt with 

CNI i the symmetry energy slope and curvature (L, K syla ) defined at p = 0.1 firT 
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I. INTRODUCTION 



The study of nuclear systems under extreme conditions is an expanding field of research involving astrophysical 
processes, laboratory experiments with rare-isotope beams, and the development of more realistic nuclear models. 
The density dependence of the symmetry energy S(p) is one of the central issues in this field 

|lH3. Although the 

symmetry energy at saturation density po is considered to be well known (J = S(po) = 33 ± 3 MeV), the different 
nuclear models present a wide range of predictions for the symmetry energy slope L = 3po[dS/dp](pa). Several 
experimental observables have been proposed to obtain a measure of L, for instance neutron-skin thickness Q, 
isovector dipolar resonances 0, Q , isoscaling in multifragmentation Q and isospin diffusion (8| . The perspective to 
obtain more stringent constraints on the value of L has motivated several studies concerning the impact of this quan tity 
on compact star structure, especially on the characterization of the core-crust transition in neutron stars [9l4ll|. In 
^v^j ' a previous work fl2j . we have addressed the model dependence of the link between L and the core-crust transition 
, properties: density pt, proton fraction Y Pit and pressure P t . It appeared that pt and Y Pt t arc unambiguously correlated 
C*~) ■ with L, while the link between L and Pt is very sensitive to the model. The transition pressure is the dominant input 
OO ' for the prediction of one of the most important crust properties: namely, the moment of inertia of the crust that may 

cn : 

affect pulsar glitches [IE [HI • ^ is then an important challenge to find a better relation between P t and laboratory 
' data. 

CNI . In this paper, we investigate the role of symmetry-energy properties other than L that could also, hopefully, be 
' related to laboratory data. To perform this study, we use as previously several Skyrmc and relativistic effective 
models, and a microscopic Brueckner-Hartree-Fock (BHF) approach. We also employ a schematic nuclear equation 
. . ■ of state (EOS), based on a density development around a reference density, that we call the generalized liquid-drop 
| model (GLDM): our objective is to establish to what extent the core-crust transition obtained with a given nuclear 
model can be reproduced with a limited set of GLDM coefficients. In other words, we investigate how the specificities 
of each model can affect the link between phenomenological constraints (i.e., the set of GLDM coefficients, which 
i could be determined from laboratory data) and the core-crust transition. 

In this paper, the core-crust transition properties are first calculated in the thermodynamic framework, with the 
transition defined as the crossing between the (3 equilibrium line and the thermodynamic spinodal border. The 
thermodynamic spinodal corresponds to the bulk instability of homogeneous nuclear matter. In fact, in a compact 
star at the core-crust transition the ground state of stellar matter changes from a clustcrized configuration (lattice 
of nuclei, or more exotic structures named pasta phases 0, [TBI . [Ill) to homogeneous matter. It has been shown 
that this transition can be very well approximated by the crossing between f3 equilibrium and the dynamic spinodal 
border [13 > which corresponds to the finite-size instability region of homogeneous matter: Coulomb and surface 
contributions make the dynamic spinodal smaller than the thermodynamic one [l8l - [20l |. Therefore, we also address 
the dynamic core-crust transition properties, which are strongly correlated with the thermodynamic ones. 

This paper is organized as follows. The generalized liquid-drop model is presented in Sec. HU We discuss the corre- 
lations existing between various GLDM coefficients and evaluate the accuracy of the core-crust transition predictions 
obtained with different GLDM expansions. In Sec. IIII1 we present in more details the analysis of the link between 
the symmetry-energy slope at saturation L and the core-crust transition, which was first presented in Ref. (l2| . In 
Sec. IIV1 we perform a correlation analysis involving coefficients others than L. This allows us to establish relations 
between selected GLDM coefficients and realistic transition properties, within a reduced range of model dispersion. 



Conclusions are presented in Sec. [V] 



II. GENERALIZED LIQUID-DROP MODEL 



Nuclear density functionals allow us to calculate the nuclear matter equation of state for any total density p = p„+p P 
and asymmetry y = (p n — p p )/(p n + p P )- They may also present a complex dependence on the density, including, 
for instance, kinetic densities, spin densities, and density gradients. Such functionals can be derived from different 
nuclear models. In this paper, we consider Skyrme and rclativistic effective models, as well as an equation of state 
based on BHF calculations. These models, presented in the Appendix, are referred to as 'complete functionals', in 
contrast with the GLDM that is presented in this section. 

The GLDM corresponds to a series expansion of the EOS around a given reference density. It is determined by 
three choices : (i) the nuclear model giving the full equation of state E(p,y), (ii) the reference density p lc f, and (hi) 
the order Af of the development. The order J\f and the reference density (p Te { ) can be chosen at will. In this paper, 
the GLDM coefficients are derived from the complete EOS E(p, y) given by various nuclear models. In a different 
context, the GLDM coefficients could be nuclear properties determined experimentally at p re f , from which we would 
extrapolate the nuclear EOS at lower or higher densities. We introduce the GLDM in order to investigate how well 
such an extrapolation allows us to reproduce the core-crust transition predicted by the complete functional; in other 
words, to what extent the role of higher-order terms can be neglected. Comparing the GLDM predictions with the 
results from the corresponding complete functional shows how much a more detailed density dependence of the EOS 
may affect these predictions. 



A. GLDM equation of state 

For a given reference density p Te { and order of development M, the GLDM energy per particle reads as: 

M x n - 

E G LDm(p, V) = / ] (ciS,n + ClV.n 2/ 2 ) T + (^kin - -E^hD > witn x = ~^ ~ ■ W 

n— 

The first term on the right-hand side of Eq.(l) contains both the kinetic and the potential contributions to the energy 
in the parabolic approximation with respect to the asymmetry y. The second term gives the contribution of the kinetic 
term beyond the parabolic approximation, as it will be explained below. We have introduced in this expression the 
GLDM coefficients cis, n an d crv.n, respectively associated with the derivatives of the energy E(p,y = 0) and of the 
symmetry energy S(p): the index 'IS' ('IV') stands for isoscalar (isovector). They are expressed as: 

ClS,n(Pref) = (3prcf)"-^-(prcf,0) (2) 

3 n S I 2 E 

ClV,n(Prcf) = (3prcf)"-^-(Prcf) with S(p) = - (p, 0) (3) 

In the case p re f = po, the lower-order coefficients are usual nuclear matter properties: Cis,o = Eq (saturation energy), 
cis, 2 = Kq (incompressibility) , Cis,3 = Qo, crv,o = J (symmetry energy), crv,i = L (symmetry-energy slope), 
c iv.2 = K sym (symmetry incompressibility), and Civ,3 = Q S ym- 

The parabolic approximation, which restricts Eq. (|TJ) to the first term, is known to be quite accurate to describe 
the EOS even at high isospin asymmetry However, it fails to reproduce the spinodal contour in the neutron-rich 
region. The reason is that the energy-density curvature in the proton-density direction must diverge at small proton 
density because of the kinetic term (see the appendix of Ref. [211]): as a result, the spinodal contour can not reach pure 
neutron matter. Instead, in the parabolic approximation, the curvature in the proton-density direction is constant, 
and leads to the unphysical prediction of unstable neutron matter. To avoid this discrepancy, we have introduced in 
Eq. (fTJ) a model-independent correction based on the non-relativistic, free Fermi gas kinetic term: 
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where r q = (3ir 2 p q ) 5 / 3 / (5n 2 ) is the kinetic density of the nucleon species q = n,p and m is the nucleon mass. As a 
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function of density and asymmetry, we have: 
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where E^ a is the parabolic part of -Ekin- In the GLDM defined by Eq. ([T]), the extra-parabolic behavior of the 
functional is sketched by the extra-parabolic behavior of Ek; n , which brings the model- independent correction Eki n — 
E^ a (second term on the r.h.s. of Eq. (JlJ). 

In the following, we will use the notation Dj^(p IC {) to identify a development of order TV around the density p re {. 
Developments of this kind are usually considered up to W=3, around the saturation density po — 0-16 fm~ 3 (see, 
e.g., Refs. [iH HU). In this work, we will also consider an extreme situation labeled E>oo- Performing an infinite 
development gives the exact value of E(p, 0) and S(p) for any density p, whatever the choice of p lc [. In practice, the 
equation of state is simply built using the exact expressions of E(p, 0) and S(p) to obtain E(p,y): 



Thus, the difference between a complete functional and its associated model is just the extra-parabolic content 
of the nuclear interaction, which has not been taken into account by the correction Eki n — E^ a defined above. Let 
us note that a parabolic approximation is often assumed in BHF calculations, to interpolate between the symmetric- 
matter and neutron-matter EOS. We will call this approach BHF para . However, by default, our BHF results include 
the correction Ek; n — Ej^ la ; in this case, the complete functional is exactly equivalent to the corresponding Doo model. 



The isoscalar and isovector GLDM coefficients obtained with all the nuclear models used in this paper are reported 
in Table |TJ These coefficients are intercorrelated, due to common constraints that the various parametrized forces 
have to verify. In particular, the fitting procedure on finite-nuclei properties, for which the average density is lower 
than po, is likely to provide effective constraints in the density region po ~0.1-0.12 fm~ 3 (this point is raised, e.g., 
in Refs. [13, HH); as a result, the coefficients defined at po for various effective models should show a tendency to 
compensate each other in order to focus the various functional predictions at lower density. To understand how a 
given coefficient (such as V) can affect the EOS properties away from the reference-density region, it is important to 
evaluate its possible connections with other coefficients. Since the equation of state of asymmetric nuclear matter is 
our main concern here, we will concentrate on the correlations between isovector properties. 

We represent in Fig. [T] different relations between GLDM isovector coefficients. On the left panel, we consider the 
coefficients taken at saturation density, crv,n(po); the different plots shows the correlation between L and </, if sym , 
and Qsym- As we have emphasized in Ref. [12j, there is a strong L-K sym correlation. The two eccentric points 
correspond to the relativistic models DDH<5I-25 and DDH<5II-30: these models include the S meson, which is generally 
associated with an atypical density dependence of the symmetry energy (see, e.g., Ref. pfij). As for the symmetry 
energy J, it tends to increase with L, although the dispersion between different models is important. The behaviour 
of the third-derivative coefficient Q aym is less universal: a clear decreasing L-Q sym correlation is obtained among the 
Skyrmc models, but the relativistic models of lower L and the BHF point arc completely out of this trend. Linear 
fits are represented on these plots; in the case of K sym and Q sym , they have been performed using only the models in 
agreement with the main trend. 

The right panel of Fig. [1] shows the relation existing between coefficients defined at two different densities, namely 
between ciy n (po) and civ, n (p=0.1fm 3 ), for n=0, 1, 2. The correlations obtained in the cases of crv,i (symmetry 
energy slope) and crv,2 (symmetry-energy curvature) reflect the shape similarities between many of the density 
functionals. The absence of correlation in the case of crv,o (symmetry energy) is due to the fact that crv,o(0-l) is 
strongly constrained in nuclear models. This discussion is illustrated in Figs. [2] (Skyrme models) and |3] (relativistic 
models). 

The left panel of each figure displays the density dependence of the symmetry energy, along with its slope and 
curvature: S(p), L(p) = (3po)[dS/dp\, and K sym (p) = {ip a ) 2 [d 2 S / dp 2 ]. Note that L(p) and K sym (p) are defined here 
with a constant factor involving po : they are not equivalent to the coefficients civ,i(/f) and crv,2(p), except at p = pa. 
We use here a constant factor in order to represent quantities proportional to the derivatives of S(p). The curves 
are shown for all models under consideration, using the complete functionals. The comparison between Figs. [5] and [3] 
shows that the relativistic models present more variability in the shape of the density functional than Skyrme models, 



E Doo (p, y) = E{p, 0) + S(p)y 2 + (E kin - E^) . 



(7) 



B. Correlations between GLDM coefficients 
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Model po Eo Kq Qq J L K synl Q By m Eoi £is,oi Koi Joi ioi -Ksym.oi 

[fm~ 3 ] [MeV] [MeV] [MeV] [MeV] [MeV] [MeV] [MeV] [MeV] [MeV] [MeV] [MeV] [MeV] [MeV] 
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TABLE I: GLDM coefficients at saturation density po and at p ro f = 0.1 fm -3 for different nuclear models. Details and 
references concerning these models are given in Appendix. BHF-1: functional fit includes calculation points in the density 
range p — [0.1; 0.35] fm~ 3 ; BHF-2: functional fit includes calculation points in the density range p = [0.05; 0.18] fm~ 3 . The 
coefficients are named following traditional notations. Isoscalar coefficients at saturation density: Eo — cis,o, Ko = cis,2, 
Qo = cis, 3- Isovector coefficients at saturation density: J = civ,o, L = civ,ii ^"sym — Civ, 2 ? Qsym = civ, 3- Isoscalar coefficients 
at p rof = 0.1 fm~ 3 : E 01 = ci S ,o, £is,oi = cis,i, Koi = cis,2- Isovector coefficients at p rci = 0.1 fm~ 3 : J i = dv,o, Un = crv,i, 

K Byul fil = Civ,2- 



as we have shown in Ref. [26| . This comes from the different ways of describing the nuclear interaction by meson 
exchange in the effective relativistic framework; in particular, the inclusion of the S meson has a strong effect on the 
density dependence of L. We can also notice that relativistic models such as GM1 and GM3, which were built for 
astrophysical purpose, have not been explicitly fitted on laboratory data. 

Considering globally the Skyrme and relativistic functionals, we can notice two common convergence regions. One 
concerns the symmetry-energy: all curves tend to cross around the density p ~ 0.11 fm -3 , taking values ^(O.ll fm~ 3 ) = 
24 ± 4 MeV. This behavior was expected, following the previous remark that finite nuclei provide fitting constraints 



at density slightly below saturation; similar observations have been made e.g. in Ref. [27]. The second convergence 
region concerns the symmetry-energy slope L(p): the different curves show a marked tendency to cross at about /?o/3. 
In contrast, no convergence effect appears for the second derivative, K sym (p). 
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FIG. 1: (Color online) Correlations between different GLDM coefficients civ.n- Left : at saturation density, relation between 
L and other coefficients civ.n(po), namely, J (top), K Byjn (center) and Q sym (bottom). Right : relation between coefficients 
civ,n defined at reference density po or p=0.1 fm -3 . Results are shown for different Skyrme models (full symbols), relativistic 
models (empty symbols), and BHF (star). BHF-f is used at saturation density, and BHF-2 is used at p=0.1 fm~ 3 . 



On the right panel of Figs. [2]and|3l we represent the density dependence of the variance between the values taken 
by different models. Denoting X(p) = {S(p), L{p), K sym (p)}, the variance is given by: 

/ = 1 N 

°x(p) = EV^)"^)] 2 ' X(P) = jjT,Xi(p) - (8) 

i i—1 

where the index i runs over the N models considered. The convergent trends arc reflected by the density dependence 
of the respective variances. In the case of Skyrme functionals, a clear minimum of as occurs at p ~ 0.13 fm -3 . With 
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FIG. 2: (Color online) Skyrme models: Left: symmetry-energy value (top), slope (center), and curvature (bottom) as a 
function of the density, calculated with complete functionals. BHF results are indicated for comparison (dots). Right: the 
corresponding variance for a calculation using the complete functionals (full line) and the GLDM up to different orders: Di(po) 
(dotted-dashed line), Z>2(po) (2dotted-dashed line), and Dz(po) (3dotted-dashed line). 
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FIG. 3: (Color online) Relativistic models: Left: symmetry-energy value (top), slope (center), and curvature (bottom) as a 
function of the density, calculated with complete functionals. BHF results are indicated for comparison (dots). Right: the 
corresponding variance for a calculation using the complete functionals (full line) and the GLDM up to different orders: Di(po) 
(dotted-dashed line), Z>2(po) (2dotted-dashed line), and Dz(po) (3dotted-dashed line). 
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rclativistic models, we observe a plateau around an inflexion point at p ~ 0.1 fm~ 3 . Note that erg is constrained to 
cancel at zero density; in such a condition, the inflexion point can be considered as a criterion of convergent trend. 
The convergence effect is even more clear in the case of the symmetry-energy slope L(p); both Skyrme and relativistic 
models present a sharp minimum of oj, at p ~ 0.06 fm~ 3 . 

These results indicate that the correlations existing between J, L, K sym , and Q S ym may be associated with effective 
constraints on the values of S(p ~ 0.12 fm~ 3 ) and L(p ~ 0.06 fm~ 3 ). To verify this, we show on the same figures the 
variances obtained when the curves X(p) are calculated using the GLDM at different orders, with reference density 
p . A minimum of its m the density range 0.11-0.13 fm~ 3 is obtained with versions Di(p ) and D 2 (po) of the GLDM. 
The version D 3 (p ) gives a clear minimum of o~l at p ~ 0.05 fm -3 . We can conclude that the correlations between 
the GLDM coefficients are related with convergence effects for S(p) and L(p) in the subsaturation density region. 

Let us remark that the convergence of the L values at p ~ 0.06 fm~ 3 is easily interpreted as a geometrical consequence 
of a constraint on S(p ~ 0.12 fm -3 ). Indeed, let us imagine a constraint fixing the values (p c , E c ) such that all models 
have to verify S(p c ) = E c . Since we also have the condition £(0) = 0, it turns out that for all models, L(p) takes 
the value L c = E c x 3po/p c at least once on the interval [0,p c ]. For instance, if p c ~ 0.12 fm -3 , E c = 25 MeV and 
po = 0.16 fm- 3 , wc have L c = 100 MeV. 

C. Comparison between GLDM and complete functionals 

We now explore the accuracy of the GLDM approximation on reproducing the core-crust transition obtained with 
the complete model on which it has been built. This procedure allows to estimate whether the core-crust transition 
can be efficiently characterized by a limited set of EOS properties determined at a fixed density. 

The discussion of the core-crust transition is performed in the thermodynamic framework, which means that, for 
simplicity, the transition is defined as the crossing between /3 equilibrium and the thermodynamic spinodal border. The 
corresponding transition density, proton fraction, and pressure are denoted, respectively, pu, Yp,tt, and Pu, where the 
index it stands for thermodynamic transition. In Sec. lIVi we will also discuss the dynamic transition defined using the 
finite-size spinodal border, in which case the index td will be used. The transition properties in both thermodynamic 
and dynamic cases are given in Table HI1 

We compare in Fig. 0] the transition properties (pu, Yp,tt, and P tt ) obtained with the complete functionals and 
different versions of their associated GLDMs. On the left panel, wc consider second- and third-order developments 
around saturation density, namely, £ ) 2(po) an d £>3(po); on the right panel, we consider the second-order development 
around a lower reference density p rc f = 0.01 fm -3 , namely, £>2(0.1). On each part, we add the results from the fully 
developped GLDM, D x , for which the only difference with the complete functional is due to extra-parabolic terms in 
the nuclear interaction. For the BHF model, the functional is equivalent by construction to its associated Doo model. 

We can see that -D2(/?o) leads generally to an important underestimation of pu, Y p ^t, and Pu- This means that, 
in a development around saturation density, terms beyond order 2 (i.e., beyond K sym ) have a large impact on the 
properties of the core-crust transition; the correlations observed between L and these properties can occur only if 
higher-order corrections are either correlated with L, or similar, for most of the functionals. The underestimation 
is strongly attenuated, but still present, with Z^po), which involves the knowledge of isovector coefficients until 
Qsym- The situation is much improved if we use a development around p rc f = 0.1 fm~ 3 . The accuracy is globally 
better with Z)2(0.1) than with D^po)] furthermore, it is nearly as good with 1)2(0. 1) as with D^, which means 
that, at p lc [ = 0.1 fm -3 , it is enough to perform the development up to second order. It is then important to relate 
experimental observables with the symmetry-energy density dependence directly in this low-density region. This 
is also appropriate, since the nucleus properties used to constrain the symmetry energy are often associated with 
subsaturation densities (neutron-skin thickness, resonances, multifragmcntation, isospin diffusion,...). 

Let us notice that effective nuclear models have an important role to play for the fine-tuning of the curvature 
properties of the EOS, crucial for the prediction of the core-crust transition. Indeed, although microscopic methods 
such as BHF are necessary to obtain reliable EOS values far away from the phenomenological constraints (in particular, 
for neutron matter and at high density), they can not be used as a reference for the EOS curvature. The numerical 
BHF results for E(p) do not allow the direct determination of second derivatives; the EOS has to be fitted (see, for 
instance, Rcf. [28]), and the resulting curvature properties are sensitive to the fitting conditions. We show in Fig. [5] 
how this affects the predictions for the core-crust transition. This figure displays the values of pu, ^,tt, and P tt 
obtained by using the BHF calculations in different ways. In the first version, BHF-1, the fit is performed on the 
density interval p = [0.1; 0.35] fm -3 : this is the version that is used to establish the GLDM coefficients at po- Ln 
order to focus on the subsaturation region, we have considered a second version, BHF-2, for which the fit is performed 
on the density interval p = [0.05; 0.18] fm -3 : this is the version we will use afterwards to define the BHF core-crust 
transition. For now, let us compare the BHF-1 and BHF-2 predictions. The GLDM coefficients at p le { = 0.1 fm -3 
have been determined for both versions, and there are significant differences in the transition properties predicted by 
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.009 


0. 


,143 


NLp6(0) 


0.069 





.014 


0.276 
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.012 


0. 


,197 



TABLE II: Density, proton fraction, and pressure at the thermodynamic transition (tt) and dynamic transition (td), for 
different nuclear models. BHF-1: functional fit includes calculation points in the density range p = [0.1; 0.35] fm -3 . BHF-2: 
functional fit includes calculation points in the density range p = [0.05; 0.18] fm -3 . 

the D2(0.1) expansion in each case. A similar contrast appears between the results obtained with the full BHF-1 and 
BHF-2 EOS. In addition to these two versions of BHF calculations, the figure also shows the complete and GLDM 
results for the Skyrme force LNS [13], the fitting procedure of which involves the BHF equation of state. From the 
span of results we obtain, it is clear that a microscopic calculation does not lead to a unique prediction for the core- 
crust transition properties; phenomenological constraints from finite nuclei will be essential to improve our knowledge 
of the low-density EOS. 
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FIG. 4: (Color online) Predictive power of the GLDM development. For different Skyrme and relativistic models, we compare 
the values of ptt (top), Y Pt tt (center) and Ptt (bottom) versus L, calculated with the complete functionals and with several 
versions of the corresponding GLDM (see text). Left: Doo, second and third order development around saturation density 
(D2,3(po))- Right: Doo and second order development around p = 0.1 fm~' ! (1)2(0.1)). 
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FIG. 5: (Color online) Comparison between exact and GLDM versions of three functionals based on BHF calculations: LNS [I?! 
(Skyrme-type force including constraints from BHF), BHF-1 (fit of the BHF EOS including densities in the range p — [0.1; 0.35] 
fm~ 3 ), and BHF-2 (fit of the BHF EOS including lower densities, in the range p = [0.05; 0.18] fm -3 ). 



III. CORRELATION BETWEEN L AND THE CORE-CRUST TRANSITION POINT 



We consider in this section the specific role of the symmetry energy slope L in the determination of the core-crust 
transition, defined here as the crossing between the line of (3 equilibrium and the thermodynamic spinodal contour. 
We summarize the study performed in Ref. [l2[, and present more details that support this previous analysis. 



A. Position of the transition point 

The most direct impact of the symmetry energy on neutron-star structure concerns the proton fraction Y p (p) in 
stellar matter, which is fixed by j3 equilibrium. For a given density, a lower symmetry energy corresponds to a lower 
proton fraction. Thus, as far as a high value of L can be correlated with a low value of the symmetry energy at 
sub-saturation, we expect higher- L models to provide lower values of Y Pt tt- This point is illustrated by Fig. El At 
fixed density p = 0.08 fm -3 , for increasing L, the /3-cquilibrium proton fraction decreases. This trend is confirmed 
and even accentuated if, instead of fixing the density, we consider the proton fraction at the transition point. The 
dispersion observed in both cases is essentially due to different values of the symmetry energy at saturation density, 



12 



0.05 r 



0.04 



o 



o 



00 0.03 - 

o 

■ 

o 

co. 

>- 0.02 



■k 



o 



o 



o 



o 



O 



0.01 



50 



100 



L (MeV) 



0.05 



0.04 



0.03- 



0.02 



0.01 




L (MeV) 



FIG. 6: (Color online) Effect of L on the proton fraction at /3 equilibrium for a fixed density p=0.08 fm 3 (left) and at 
thermodynamic spinodal crossing (right), for different nuclear models: Skyrme (full symbols), relativistic (empty symbols), 
BHF-2 (star) and BHF-2 para (asterisk). 



as will be discussed in the following. 

Let us now consider the impact of L on the transition density p u , illustrated by Fig.[7J The correlation between L 
and pu is a well-known feature 0; however, its explanation is less intuitive than in the case of the L-Y pM correlation. 
Furthermore, it can not be explained just as a consequence of the behavior of Y p ttl as we see on the left panel of 
the figure: even for a fixed proton fraction Y p = 0.025, the spinodal border shows a clear decreasing correlation with 
L. This feature can be understood as a consequence of the strong link existing between L and the energy-density 
curvature of neutron matter taken at symmetric spinodal density p s : 

2 T 1 V- x s~ 2 
3po 3p ^ (n - 2)! 

with a; s = (p s — po)/(3po)- Since all models yield a symmetric spinodal density close to 0.1 fm -3 , so that x s ~ —1/9, 
the term n = 2 is nearly canceled: this reinforces the dominance of L in the determination of Cnm.s- 

Figure [5] illustrates how L affects independently the proton fraction at /3 equilibrium and the spinodal contour in 
the neutron-rich region. Due to the typical geometry of these respective lines, the two effects reinforce each-other, 
leading to a robust correlation between L, pu and Yp,tt- 
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(9) 



B. Pressure at the transition point 



The link between L and the core-crust transition pressure Pu is more problematic. In order to make this link 
explicit, let us write the pressure in the GLDM framework: 
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FIG. 7: (Color online) Effect of L on the density at thermodynamic spinodal border for a fixed proton fraction Y p — 0.025 
(left) and at /3 equilibrium (center). Right: correspondence between L and the energy-density curvature of neutron matter at 
symmetric spinodal density, Cnm,s (see text). Results are shown for different nuclear models: Skyrme (full symbols), relativistic 
(empty symbols), BHF-2 (star) and BHF-2 pa ra (asterisk) (both BHF results are identical for neutron matter). 




FIG. 8: (Color online) Simultaneous effect of L on the proton fraction at ft equilibrium and on the thermodynamic spinodal 
border of neutron-rich matter. Left: spinodal contours and /3 equilibrium for three Skyrme models: BSkl7 (L=36 MeV), SkI6 
(L=60 MeV), Ro- (L=86 MeV). Right: correlation between ptt and Y v ,tt for different nuclear models: Skyrme (full symbols), 
relativistic (empty symbols), BHF-2 (star) and BHF-2 para (asterisk). 



From this expression, we expect that for a given density the pressure of neutron-rich matter should increase with L, 
which is the leading coefficient. This trend appears on the left panel of Fig. [9], representing the relation between L and 
the pressure of pure neutron matter, Pnm, at p=0.08 fm -3 . Thus, a positive correlation between L and Ptt should 
be obtained if we could neglect the density shift due to the L-p tt correlation, as well as the effect of higher-order 
coefficients. However, as it can be seen on the right panel of the figure, the results for P tt (L) present an important 
dispersion and we cannot extract a clear correlation, although a decreasing trend can be observed among Skyrme 
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FIG. 9: (Color online) Impact of L on the pressure of neutron-rich matter at sub-saturation density, for different nuclear 
models: Skyrme (full symbols), relativistic (empty symbols), BHF-2 (star) and BHF-2 para (asterisk). Both BHF-2 results are 
identical for neutron matter. Left: pressure of pure neutron matter, for a fixed density typical of the transition: p=0.08 fm -3 . 
Right: thermodynamic transition pressure. 



models. Four eccentric points close to L=60 MeV weaken this correlation between Skyrme models: they correspond 
to atypical relations between L and K sym , which also affect the plot Pnm(0.08 fm -3 ). 

The lack of correlation between L and Pu when independent models are considered results from a delicate balance 
between opposite effects, as we have discussed in Ref. [l2j]. This is shown by separating the different contributions 
we can estimate from the GLDM formula. We distinguish two kinds of contributions to the variation dP u /dL: (i) 
variations occurring at a fixed density (p, y) , resulting only from the modifications of the coefficients in Eq. (|10l) , 
which defines P(p 7 y); and (ii) variations due to a shift (Sp 7 Sy), for a fixed expression of P(p,y), i.e. frozen values of 
the coefficients in Eq. (|10[) . The contributions of the first kind come from the explicit L dependence of Eq. (|TU)) . and 
from correlations between L and higher-order coefficients civ.n- In practice, we will consider the following terms: 
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where the index e means that the modification concerns the expression of the pressure, and the number gives the order 
of the modified coefficient. The contributions of the second kind, resulting from the density position of the transition 
point, are characterized by the index p. We distinguish the respective effects of total density and asymmetry: 
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The variations of the quantities depending on L arc fixed empirically, using as a reference the correlations that are 
observed between different models. From linear fits, we extract Spu/SL=-3.84: x 10 -4 MeV -1 fm -3 , Syu/SL=6.08 x 
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FIG. 10: (Color online) Estimation of different contributions to the variation of P t t with L (see text). Left: separated 
contributions, due to the transition position shift (indices pi, p2) and to the variation of the GLDM coefficients crv.n in the 
expression of P (indices el, e2, and e3 for n=l, 2, and 3 respectively). Right: sum of the contributions, considering the 
contribution of 5civ,„/5L up to order n=l, 2 and 3 (respective index: totl, tot2, tot3). The horizontal lines indicate the region 
where the pressure variation is compatible with zero within the estimated uncertainty. Full symbols: Skyrme models; empty 
symbols: relativistic models. 



10~ 4 MeV -1 , SK aym /SL=3.33, 5Q s /SL=-6.63. Note that, in the case of Q s , the correlation with L is observed only 
within the Skyrme models, which are used to perform the linear fit. 

These various contributions are represented on Fig. 1101 It appears that the contribution of the asymmetry shift is 
quite marginal, and the overall SPu/SL essentially results from the balance between three terms: [5P t t/ 8L] cl , which 
is large and positive, is compensated by the conjugated effect of [5Pu/SL] o2 and [SPtt/5L] pl . These two negative 
contributions, e2 and pi, are of the same order of magnitude: this means that the correlation L-K sym has the same 
importance as the correlation L-p u in explaining why we do not observe an increasing correlation P t t{L). 

In addition, we can see that the term (5 Ptt) / '5L) e ^ due to the L-Q sym relation brings an additional negative 
contribution, but of much lower magnitude than the term e2. This result means that, although the third-order term 
of the GLDM has a strong impact on the absolute value of Ptt , as it was observed on Fig. 01 the L-Q sym correlation 
is not crucial in the determination of SPu/SL; in other words, the third-order correction does not depend strongly 
on L. On the other hand, let us notice that the strong dispersion of Qsym values in the case of relativistic models is 
bound to cause a strong dispersion in P t t{L). 

To summarize, if we characterize the L-P u relation using a GLDM development around saturation density, we can 
identify three effects that are crucial for the determination of SP/SL: (i) the explicit L dependence of the pressure 
given by Eq. (fTUj) ; (ii) the L-p u correlation and (iii) the L-K sym correlation. These different contributions compensate 
each other. For some models (those of higher L), the GLDM predicts a decreasing P tt (L); for others (those of lower 
L), an increase would be obtained. It is interesting to note that SP/SL cancels in the interval of the most realistic L 
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FIG. 11: (Color online) Typical variation for the values of J, K ayul (L) and Kq. We represent these quantities for different 
Skyrme models (full symbols), relativistic models (empty symbols), and BHF-1 (star). The lines show the different values 
which are used for the standard schematic model (see text). 

values, namely 50-80 MeV. By estimating an uncertainty of about 20% on the slopes of the L-p u and L-K sym linear 
fits, we obtain an error bar of ±0.02 fm -3 on 5P/5L, which appears compatible with zero throughout this interval. 
These results are not a quantitative prediction on the evolution of the transition pressure with L\ however, they show 
that the link between L and Ptt cannot be deduced from qualitative arguments, and therefore it is not soundly based. 
The relation between L and Pu is in fact very sensitive to model-dependence, as it will be further discussed in the 
following. 

C. Predictions of a standard GLDM 

In order to study how the different GLDM coefficients can affect the core-crust transition, we will make use of a 
schematic model corresponding to a D^{po) expansion with typical values for the different coefficients. The choice 
of these values is illustrated on Fig. 111! which gives a graphical representation of the saturation properties of the 
different nuclear models considered in this paper. The lines indicate the intervals of coefficients attributed to the 
typical GLDM that we are now constructing. We define a reference D^(po) model characterized by the following 
parameters: 

po = 0.16 fm" 3 ; K = 240 MeV ; Q = -350 MeV ; J = 31 MeV . 

The coefficient L varies in the interval [40; 100] MeV, and determines K sym and Q sy m according to the relations: 

K sym (L) = a K xL + b K ; (a K ; &k) = (3.33; -281) 
Q S (L) = a q xL + 6 Q ; (a Q ; b Q ) = (-6.63; 765) 

The predictions of this standard model for p tt (L), Y pM (L), and Ptt(L) are represented in Fig. [T2] As expected, we 
obtain a clear decrease of p u and Y p ^t with L, while the evolution of Ptt is quite flat. In the following, we will observe 
how these curves evolve when some of the standard EOS properties arc modified within a realistic interval. We will 
modify separately the symmetry energy J, the incompressibility Kq, and the L-K sym relation. On the left panel of 
Fig. [T2J we show the effect of varying J between 29 and 33 MeV. On the right panel, we show the effect of varying 
Kq between 220 and 260 MeV. On the central panel, we modify the linear relation between L and K sym by adopting 
different values of (ok! &k) : a softer version (ok; friOsoft = (2-8; —235), and a stiffer one (ax; &K)stift = (3-8; —320). 

The main features to be noted concerning the results of the standard schematic model and how they are affected by 
typical variations of the GLDM coefficients are the following, (i) The qualitative behavior of pu and Y p ^t is maintained: 
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FIG. 12: (Color online) Predictions of the standard schematic model for the transition density, proton fraction and pressure, 
compared with complete effective models: Skyrme (full symbols) and relativistic (empty symbols) models. Left: varying the 
symmetry energy at saturation, J. Center: varying the relation K aynl (L). Right: varying the incompressibility at saturation 
K . 



although the absolute value can be affected by different aspects of the functional, they always unambiguously decrease 
with increasing L. (ii) The qualitative behavior of Ptt is very sensitive to the values of the GLDM coefficients: the 
application of a very moderate variation, inside a realistic model uncertainty, leads to opposite predictions: P tt (L) 
either increases or decreases, and most often it is quite flat. These two conclusions confirm the previous analysis. 
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FIG. 13: Relation between the dispersion of Y P: tt{L) and the dispersion of GLDM properties (J, Ko, AK ayal (L)) associated 
with the different models (see text). Nuclear models: Skyrme (full symbols), relativistic (empty symbols), and BHF-2 (star). 



IV. PARAMETRIZATION OF THE CORE-CRUST TRANSITION 



In this last section, we explore the possibility to reduce the model dispersion in the prediction of the core-crust 
transition, by taking into account the effect of coefficients others than L. First, we check to what extent the dispersion 
in the L dependence of ptt, Y PtU , and Ptt can be attributed to specific GLDM coefficients, such as the symmetry energy 
at saturation J, the incompressibility Kq, or the quantity AK sym (L) = A' sym — (ax x L + b^). This last quantity 
characterizes the eccentricity of the model with respect to the typical relation K sym (L) = (ax x L + b^). On the 
basis of these results, we propose to fit the core-crust transition properties by a linear dependence on pairs of GLDM 
coefficients. This idea is applied to the thermodynamical transition, which has been the framework of our analysis, 
and to the dynamic transition, which is the best approximation to the realistic core-crust transition. 



A. Role of the first GLDM coefficients in the dispersion 



As we have seen in the previous section, the correlations between L and the core-crust transition properties ptt, 
Yp^t, and Ptt suffer from a certain amount of dispersion when different kinds of models are considered. This effect 
is particularly harmful for the prediction of the transition pressure; with the link between L and Ptt being very 
sensitive to the details of the functional, the model dispersion destroys the possibility to deduce the value of Ptt from 



19 




FIG. 14: Relation between the dispersion of ptt{L) and the dispersion of GLDM properties (J, Ko, AK syln (L)) associated 
with the different models (see text). Nuclear models: Skyrme (full symbols), relativistic (empty symbols), and BHF-2 (star). 

a measurement of L. In order to look for model properties that may be responsible for this dispersion, we first use 
the following procedure: 

• We call M the GLDM model property whose effects are investigated: namely, J, Kq or AK sym . 

• We represent the quantities X(L) calculated with the different models, with X=Y p tt (Fig. [T3\i , p tt (Fig. HPi]) . 
and P u (Fig. H5]). 

For each M property (J, Kq or AK sym ), we check whether the diagram X(L) can be separated in two regions 
associated with larger/smaller values of M. For this: 

• we use a trial frontier, namely a straight line A: X\ = a x L + b 

• we calculate the distance of each data point i to this frontier: [AIJa^ = I, - (a x Li + b). 

• we vary the frontier A in order to obtain the best correlation for the diagram [AX]\(M). 

In this way, a different line A is defined for each model property M under study. This is clearly seen in Figs. 1131 
[T4l[T5l where the frontier line in the three top graphs, associated, respectively, with M = J, AK sym , Ko, differs from 
one graph to the other. 
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FIG. 15: Relation between the dispersion of Ptt(L) and the dispersion of GLDM properties (J, Ko, AK sym (L)) associated 
with the different models (see text). Nuclear models: Skyrme (full symbols), relativistic (empty symbols), and BHF-2 (star). 



If the dispersion of the data points can be mainly attributed to differences in their respective M values, the diagram 
[AX]\(M) must present a clear correlation. In this case, it would be sufficient to know the values of L and M of a 
given model to predict accurately the corresponding X value. On the contrary, if the diagram [AX]a(M) appears 
completely uncorrelated, we can deduce that the coefficient M under study is not responsible for the dispersion of 
the data Xi(Lj), In several cases, we find an intermediate situation, in which the correlation of [AX]a(M) is weak 
but allows to associate eccentric values of M with large values of AX. 

The most favorable situation appears in Fig. Q21 with the effect of J on AY Pytt . This effect was expected : Y pM 
depends on the value of the symmetry energy at sub-saturation, which is well correlated with L as long as the different 
models have a similar symmetry energy at po- J also appears to affect the values of p u (Fig. Q3]) and P t t (Fig. IT5"]) . 
although the correlation [AX]a(J) is weaker in these two cases. Let us now consider the eccentricity of the K sym 
behavior, namely, AK sym . It has a clear effect on the relation P U (L), as we see in Fig. [15] this confirms the analysis 
of the previous section, where we have underlined the role of the L-K sym correlation in the link between L and Pu- 
On the other hand, AK sym is uncorrelated with the position of the transition, p u and Y Ptt t- Finally, the isoscalar 
incomprcssibility Kq has no clear effect on the core-crust transition. It appears completely uncorrelated with the 
values of Y p ^t and Pu- A weak correlation appears with p u , but this does not affect significantly the quality of the 
L-pu correlation. In the following, we will concentrate exclusively on the role of isovector coefficients. 
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FIG. 16: Comparison between the thermodynamic (dashed line) and dynamic (full line) spinodals. The dotted line repre- 
sents the /^-equilibrium EOS and the square and dot define the crust-core transition within, respectively, the dynamical and 
thermodynamical spinodal. 
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FIG. 17: (Color online) Comparison between the transition taken at dynamic (id) and thermodynamic (tt) spinodal for different 
Skyrme (full symbols) and relativistic (empty symbols) models: density (left), proton fraction (center), and pressure (right). 



B. Prediction of the dynamical core-crust transition 

Until now, we have studied the quantities X tt = {ptt, Yp,tt> Ptt}, which are defined by the crossing between the j3- 
cquilibrium condition and the thermodynamic spinodal. This framework allowed us to emphasize the analytical role of 
bulk GLDM coefficients in the transition. However, realistic descriptions of the core-crust transition involve stability 
comparison between homogeneous matter and clusterized matter. Equilibrium calculations have been performed, e.g., 
in Refs. 0, Il5l - [r7| : it has been verified that the resulting transition can be very well approximated by the crossing 
between the /3-cquilibrium condition and the dynamic spinodal. In Fig. 1161 we illustrate the differences between the 
thermodynamic and dynamic spinodals and identify with a square (dynamic) and dot (thermodynamic) the crust-core 
transition, defined by the crossing of the /3-equilibrium EOS and the spinodal. We denote X t d = {ptd, Yp,td, Ptd} as 
the quantities taken at this dynamic spinodal border. We have obtained these quantities within the effective Skyrme 
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and rclativistic approaches. The dynamic spinodal has not been calculated in the BHF framework, since the BHF- 
based density functional does not include the density-gradient terms needed to modelize the surface effects associated 
with finite-size density fluctuations. 

To check that the thermodynamic framework effectively reflects the correlations between the GLDM coefficients 
and the core-crust transition, we have to make sure that the transformation from X u to X tc i does not destroy these 
correlations. This is verified in Fig. 1171 where we plot the dynamic results as a function of the thermodynamic ones: 
we observe that these quantities are strongly correlated. In the following, we will extract some empirical relations 
between the GLDM coefficients and the dynamic core-crust properties, {ptd, Y p jd, Ptd}- 

A systematic analysis of the effect of the isovector coefficients J, L, and K sym on the transition properties is 
presented in Table Hill Our previous study has shown that the relation X a {L) is affected by atypical values of J 
and -Ksym associated with the various models. To explore this effect, we have performed two-dimensional fits of the 
transition data: 

X ti (M 1 ,M 2 ) = a 1 M 1 +a 2 M 2 + b, (16) 

where Xu = {ptt,Y p jt, Ptt] ptd,Y Pj td, Ptd} (dynamic or thermodynamic transition), and Mi are two of the isovector 
GLDM coefficients: J re f = S(p re f), L Ie f = 3p Ic f[dS/dp](p Te {), and A sym , rcf = (3p TC f) 2 [d 2 S/dp 2 ](p TC f). We have 
considered the saturation coefficients J, i, and A sym , as well as coefficients at the reference density p Te { = 0.1 fm -3 , 
denoted Joi, Lqi, and if S ym,oi- Table Hill gives the root mean square (rms) of residuals associated with the different 
fits, indicating the relevance of the respective combinations of coefficients in the determination of Xu- 

It appears that p t i is well correlated with L; no significant improvement can be obtained by considering pairs 
of coefficients. The quality of the L-p t i correlation can be understood as a consequence of Eq. ([5]), as discussed in 
Section flH Al In the cases of Y Pt u and Pu however, the predictions can be considerably improved by using combinations 
of coefficients. As expected, Y p _u is very well correlated with a combination of J and L, and this result is still improved 
using a combination of Joi and Lqi- On the other hand, the values of Tab. IIIII indicate that combinations of L and 
-Ksym are not relevant to determine Y Pj u- The transition pressure, instead, presents improved correlations with two 
kinds of parameter combinations: either J and L at saturation density, or Lqi and K'sym.oi; the latter leads to the 
smallest rms of residuals. 

We show in Fig.[TS]the most significant correlations obtained between the dynamic transition properties and different 
combinations of the GLDM coefficients. The quantities p t d, Y p ^d, and Ptd are shown as a function of L and as a 
function of selected coefficient combinations. For the transition density, we verify the good L-p t d correlation; we also 
consider two coefficient combinations, [J — 0.558 x L] and [Loi + 0.426 x A" sym .oi], which lead to a similar dispersion, 
lower than ±0.004 fm" 3 . The situation is different for the proton fraction Y p ^d and the pressure Ptd, for which 
the correlation with L is not so good. The selected coefficient combinations lead to clearly improved correlations. 
Concerning the proton fraction, the combinations [J — 0.172 x L] and especially [Joi — 0.160 x Joi] present excellent 
correlations with Y Ptt d- As for the transition pressure, the combinations [J — 0.127 x L] and [Joi — 0.343 x -K S y m ,oi] 
allow to considerably reduce the data dispersion, and show unambiguous correlations with P t d', this is especially true 
for the second combination (coefficients extracted at 0.1 fm -3 ), for which the typical model dispersion for P t d becomes 
± 0.033 MeV.fm 3 instead of ± 0.085 MeV.fm 3 in the case of P td (L). The linear fits represented in Fig. 1X51 are: 



p td (L) = (-3.75 x 10" 4 x L + 0.0963) fm -3 (17) 

p td (J, L) = (7.46 x 10" 4 x [J - 0.558 x L] + 0.0754) fm -3 (18) 

p td (£oi,-K 8 ym,oi) = (3.23 x 10- 4 x [Loi + 0.426 x X sy m,oi] + 0.0802) fm" 3 (19) 

Y p j d (L) = -3-15 x 10 -4 x L + 0.0461 (20) 

Y pM {J, L) = 2.69 x 10" 3 x [J - 0.172 x L] - 0.0290 (21) 

Y pM (Joi,L i) = 3.30 x 10" 3 x [J 01 - 0.160 x L 01 ] - 0.0253 (22) 

P td {L) = (-2.42 x 10~ 3 x L + 0.465) MeV fm" 3 (23) 

P td (J, L) = (3.36 x 10" 2 x [J - 0.127 x L] - 0.474) MeV fm" 3 (24) 

Ptd(L Q i,K sym . Q1 ) = (9.59 x 10" 3 x [L 01 - 0.343 x A; ym . 01 ] - 0.328) MeV fm" 3 (25) 



The relations involving GLDM coefficients at p=0.1 fm" 3 allow to predict the core-crust transition density, asymmetry 
and pressure within a reasonable model uncertainty. This indicates that exploring the sub-saturation properties 
of nuclei in order to constrain directly these low-density coefficients could allow to modelize core-crust transition 
properties that do not depend strongly on the type of nuclear functional that is used. 
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TABLE III: Root mean square of residuals associated with two-variable linear fits X« = ai X Mi + a 2 x M 2 + b, where 
Xu = {ptt, Y p ,u, Ptt; ptd, Y p ,td, Ptd} and Mi are two of the isovector GLDM coefficients: J rc f = S(p re f), L rc c = 3p rc {[dS/dp}(p re f) 
and -Ksym.rcf = (3p re i) 2 [d 2 S/dp 2 ](p rc f). J, L and -ft" sym are taken at p Te t = po', Joi, Loi and 7f S ym,oi are taken at p Te i = 0.1 
fm -3 . The results for single- variable linear fits in function of L and Loi are also shown. 

V. CONCLUSION 

In this paper, we have studied to what extent the core-crust transition properties can be predicted by using a 
reduced set of phenomenological constraints. 

Different kinds of nuclear models have been compared: Skyrme and relativistic effective models, and a microscopic 
BHF approach. It is seen that the EOS obtained in the relativistic cases present much more variability in their density 
dependence than the Skyrme ones, due to different ways of describing the interaction; notably, a softer symmetry 
energy is obtained with density-dependent couplings, and the inclusion of the delta meson leads to an atypic density 
evolution of the symmetry energy slope, as was noticed in Ref. ,2fi|. A more regular behavior is observed within the 
Skyrme sets of parameters; this has the drawback to bring possibly spurious correlations between the EOS properties at 
saturation and subsaturation density. As for the BHF calculations, although microscopic approaches are indispensable 
to provide realistic predictions for the EOS away from the phenomenological constraints, their predictions for the 
core-crust transition properties are very sensitive to the numerical fit of the EOS, which is necessary to determine its 
curvature. However, other ways could be followed to constrain the curvature properties from microscopic approaches, 
such as the study of Landau parameters [2!| [3(| ■ 

To relate the predictions of the different models with their properties at a fixed density, we have introduced 
a generalized liquid-drop model (GLDM) which consists of a density development of the EOS around a reference 
density p re f, up to a chosen order. When p re f is the saturation density, we have seen that a development up to 
order 3 is necessary to get reasonably close to the thermodynamic core-crust transition properties predicted by the 
complete functionals. This means that the correlations that can be observed between the symmetry-energy slope 
at saturation L and the core-crust transition properties are subjected to further correlations existing between the 
various coefficients of the GLDM. Such correlations reflect two kinds of effects: (i) the possible existence of effective 
constraints at subsaturation densities, arising from the nuclear data used in the fit of most of the effective models, and 
(ii) the regularity of the functional shapes, which depend on the construction of the model. We have also considered 
a development at a lower reference density, p rc [ = 0.1 fm . This approach has the advantage of reducing the 
model dependence arising from specific functional shapes, and of focusing on a density region closer to most nucleus 
observations. A development at order 2 around p — 0.1 fm -3 allows us to characterize the thermodynamic core-crust 
transition within nearly the smallest uncertainty allowed by a GLDM approach. This smallest uncertainty is given 
by the infinite development D^, which gives the best GLDM approximation of the complete functional by neglecting 
only the extra-parabolic terms in the isospin dependence of the nuclear interaction. 

We have also presented a more detailed study of the relation between L and the core-crust transition, the conclu- 
sions of which confirm our previous analysis [12j. Namely, the core-crust transition density pu and proton fraction 
Y Pt tt appear clearly correlated with L, despite the variety of models, while the link between L and the transition 
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FIG. 18: (Color online) Correlations between the dynamic core-crust transition properties and GLDM coefficients (see text). 
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pressure Pu is much more sensitive to model dependence. Indeed, the impact of L on the transition pressure involves 
several opposite contributions, which tend to compensate each other; thus, it is not possible to establish a qualitative 
prediction for the evolution PttjL\ which can change sign depending on the model (see, for instance, the opposite 
predictions presented in Refs. [ If! [ilj ) . 

To explore the possibility to overcome the model dispersion and predict the core-crust transition properties from a 
reduced set of nuclear constraints, we have searched those GLDM coefficients others than L that play a major role 
in the determination of this transition: J and -Ks ym were found to have a significant responsibility in the observed 
dispersion. Finally, we have addressed the case of the dynamic core-crust transition, given by the crossing between 
the dynamic spinodal and the (3 equilibrium. This corresponds to a realistic approximation of the actual core-crust 
transition, and takes place at lower density than in the thermodynamic approach. We have verified that the dynamic 
transition is related to the GLDM coefficients by similar correlations. The L-p t( i correlation is quite good, and cannot 
be significantly improved by considering other coefficients; however, the predictivity of Y Ptt d and Ptd is considerably 
better in terms of selected pairs of coefficients. An excellent correlation appears between Y Pit d and a combination 
of J and L, and it is even better using a combination of Joi and Loi (coefficients defined at p Te t = 0.1 fm~ 3 ). 
Furthermore, the model dependence in the prediction of the transition pressure can be considerably reduced if we 
consider a combination of Loi and -Ksym.oi- In this case, an unambiguous correlation is obtained within all the variety 
of Skyrmc and relativistic models considered. 

To conclude, it appears that an accurate determination of the first three GLDM coefficients at p lc f = 0.1 fm -3 would 
allow a prediction of the core-crust transition properties that do not depend much on the model construction. This 
gives a strong motivation to focus on the relation between nuclear observables and GLDM coefficients at subsaturation 
density. It will become possible to use phcnomenological nuclear models to restrict the range of the core-crust transition 
properties in neutron stars; this would have an impact on the interpretation of astrophysical observations, and on the 
possible scenarios to explain phenomena such as pulsar glitches fl3l [l4| . As discussed in (l3j , the transition pressure 
is an essential input to infer the neutron-star mass-radius relation from glitch observations. Since the mass-radius 
relation predicted by a given EOS is mainly determined by its high-density region, an accurate prediction of the 
transition pressure would also constrain the high-density EOS. 



Appendix A: Nuclear models 



In this appendix, we give an overview of the different types of nuclear models whose properties are compared in the 
present work. Each of these models gives at every density the energy of symmetric matter and the symmetry energy 
(in contrast with the GLDM exposed in the Section Hi Al which is based on a density expansion). 



1. Skyrme-like effective models 

The local Skyrme interaction [3l[ allows us to define an energy density functional H(r) such that the total energy 
for a system of nucleons in a Slater determinant | ip > reads as 

(ip\H\ip) = J H(r)d 3 r . (Al) 

For homogeneous, spin-saturated matter with no Coulomb interaction, the Skyrme energy density functional is 
composed of four terms: 

H = JC + H +H 3 + n cS (A2) 

In this expression, K, is the kinetic-energy term, 7fo is a density-independent two-body term, TI3 is a density-dependent 
term, and TL e g is a momentum-dependent term [32j: 

K = f- , (A3) 
2m 

n = C oP 2 + D pl , (A4) 

n 3 = Csp^ 2 + D s p°pl , (A5) 

Hoft = Coffpr + D cS p 3 r 3 . (A6) 

We have introduced the isoscalar and isovector particle densities, p and p 3 , as well as kinetic densities, r and t 3 : 



P = Pn + P P ; i 
P3 = P 



(A7) 
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where, denoting q the third component of the isospin (n for neutrons and p for protons), the kinetic densities are 
defined by r q = {k 2 ) q . The coefficients C and D, associated respectively with the isoscalar and isovector contributions, 
are linear combinations of the traditional Skyrmc parameters: 

Co = 3^o/8 , 

D = -t (2x + l)/8 , 

C 3 = t 3 /16, , 

D 3 = -t 3 (2x 3 + l)/48 , [A *> 

Ceff = [3ii+t 2 (4a:2 + 5)]/16 ) 

Doff = [ta(2a;3 + l)-«i(2xi + l)]/16. 

In this paper, we have considered 21 Skyrmc paramctrizations, commonly used in the literature, chosen in order to 
cover a wide range of L values while presenting acceptable saturation properties. Traditionally, Skyrmc parameters 
are fitted in order to reproduce selected nuclear properties measured in a set of nuclei: basically masses and radii, plus 
several other input with increasing level of sophistication. SV [33| is among the earlier parametrizations. SGII [IH], for 
which s pin p roperties have also been used as constraints, can reproduce isos pin effects in giant dipole resonances. R CT 
and G CT [35( consider spin-orbit splitting in 16 and surface widths. SkMP [36[ was built to improve the fit of 208 Pb 
charge distribution. The series SkI2, SkI3, SkI4, SkI5 [37| and SkI6 [38| include constraints on the isotope shifts of the 
charge radius in Pb and Ca. SkO [39( further considers isotopic evolution of two-neutron separation energies in Pb. In 
addition to nuclear data constraints, many Skyrme forces include in their fitting procedure the neutron matter EOS 
from microscopic calculations: the objective is to obtain a reliable behavior of the density functional at high isospin 
asymmetry, especially for astrophysical applications. RATP [40| was the first parametrization using this procedure, 
including the neutron matter calculation by Friedman and Pandharipande [4l[ . The Skyrme- Lyon forces SLy230a, 
SLy230b [H| , SLv4 _|42| and SLylO [H] use the pure neutron matter equation of state UV14+UVII by R.B. Wiringa 
ct al (44[. NRAPR [45[ (Non-Relativistic APR) stands for the Skyrme interaction parameters obtained from a fit of 
the APR equation of state (Akmal-Pandharipande-Ravcnhall, Ref. [Hj]). LNS (47j is based on Brueckner-Hartrce- 
Fock calculations of infinite nuclear matter at different values of isospin asymmetry. The Bruxcllcs-Skyrme forces 
BSkI4 [48| . BSkl6 (4!| and BSkl7 [5(| include the Friedman and Pandharipande calculation of neutron matter [4lj . 
and a HFB treatment of pairing effects in order to improve mass predictions in the neutron-drip region. 



2. Relativistic effective models 



In this paper, we consider two kinds of relativistic effective approaches: RMF models, which have constant coupling 
parameters described by the Lagrangian density of non- linear Walecka models (NLWM), and DDH models with 
density-dependent coupling parameters. In each case, we consider models including or not the 8 meson, which have 
been introduced to include in the isovector channel the same symmetry existing already in the isoscalar channel 
with the meson pair (a, oS) responsible for saturation in RMF models |51j . The presence of the 8 meson softens the 
symmetry energy at subsaturation densities and hardens it above saturation density. The RMF parametrizations we 
use are NL3 [Hf, TM1 GM1, GM3 H, FSU [H, NLwp [H, and NLpS H3- The DDH parametrizations are 
TW [H, DD-ME1, DD-ME2 [H, and DDH<5 0. The models NLp<5(0) and NLp£(2.5) introduced in H3 have the 
same isoscalar properties and the same symmetry energy at saturation; however, the last model includes the S meson 
with gs = 2.5, while in the first one, the 6 coupling was set to zero. We have also introduced the parametrization 
NLp(5(1.7) with a weaker <5-meson coupling (1.7 instead of 2.5). The model NLo;p(025) includes a ujp non-linear term 
in the Lagrangian as in (5|| with strength A v = 0.025. The parametrization DDHi5I-25 introduced in [6(J has a quite 
low symmetry energy at saturation (25 MeV); therefore, we also consider the parametrization DDH5II-30 where the 
p-meson coupling was adjusted so that, at saturation, the symmetry energy is 30 MeV, and all the isoscalar properties 
are kept fixed. 

The relativistic approach is based on a Lagrangian density given by: 



The nucleon Lagrangians read as 
with 



i— p,n 

iD» = idV-TvVf —t • b^ 1 , 

M* = m-T s (f>-T S T -5 , 



(A9) 

(A10) 

(All) 
(A12) 
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where t is the isospin operator. We use the vector symbol to designate a vector in isospin space. 

The isoscalar part is associated with the scalar sigma (a) field 4> and the vector omega (u) field Vn, while the isospin 
dependence comes from the isovector-scalar delta (8) field S 1 and the isovector-vector rho (p) field 6^ (where p, is a 
space-time index and z an isospin-direction index). The associated Lagrangians are: 



^ = + l<v^ + y g 4 v (v,vn 2 , 

where f2 M „ = d^Vv — d v V^, B^ u = d^b v — d v b p — T p (6 /1 x b v ), and Tj and rrij are, respectively, the coupling parameters 
of the mesons j = s, v, 8, p with the nucleons and their masses. The self-interacting terms for the a meson are included 
only for the NL3 and NL5 parametrizations, with k and A denoting the corresponding coupling constants. 

The density-dependent coupling parameters T s , T v , and T p , are adjusted in order to reproduce some of the nuclear 
matter bulk properties, using the following paramctrization: 



with 



where x = p/p sa t and 



Ti(p) = Ti(p sat )fi(x) , i = s, v ( A1 3) 



1 + Ci{x + diY 



T p {p) = T p (p sat ) cxp[-a p (x - 1)] . (A15) 



The values of the parameters m,, Tj, a,i, bi, Ci, and di, i = s,v, p, for TW and DD-ME2 are, respectively, given in [58[ 
and [6l[ and for DDH<5 in [6(J[62j]. In this last case, the parametrization for the 5- and p-coupling parameters is also 
given by (|A13j) with 

fi(x) = aiexp[-bi(x - 1)] - Ci(x - di) , i = p, 5. 
The Ti coupling parameters are replaced by the gi coupling constants in the NL3 and NL<5 models. 



3. The BHF approach of asymmetric nuclear matter 

The BHF approach of asymmetric nuclear matter [63, 64] starts with the construction of all the G matrices describing 
the effective interaction between two nucleons in the presence of a surrounding medium. They are obtained by solving 
the well-known Bethe-Goldstone equation 

^-*TlT2;T3T4 (^0 — ^ / TlT2,'T3T4 ~l~ / , ^Tl T2 ; 7~i T 7 ' \ '. G T . Tlj - T ^ T4 (tj) , (A16) 

— ' U! — Ci — 6j + If) 

ij J 

where r = n,p indicates the isospin projection of the two nucleons in the initial, intermediate, and final states, V 
denotes the bare NN interaction, Q TiT - is the Pauli operator that allows only intermediate states compatible with 
the Pauli principle, and w, the so-called starting energy, corresponds to the sum of non-relativistic energies of the 
interacting nucleons. The single-particle energy e T of a nucleon with momentum k is given by 

e T (k) = ^-+-RB[U T (k)] , (A17) 
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where the single-particle potential U T {k) represents the mean field "felt" by a nucleon due to its interaction with the 
other nucleons of the medium. In the BHF approximation, U (k) is calculated through the "on-shell energy" G-matrix, 
and is given by 

U T (k) = J2 (**' I G T t>;tt>(u = e T (k) + e T ,(k')) | kk') A (A18) 

T ' \k'\<kp T , 

where the sum runs over all neutron and proton occupied states and where the matrix elements are properly antisym- 
metrized. We note here that the so-called continuous prescription has been adopted for the single-particle potential 
when solving the Bethe-Goldstone equation. As shown in Rcfs. (65l. |66|. the contribution to the energy per particle 
from three-hole line diagrams is minimized in this prescription. Once a self-consistent solution of Eqs. (|A16|) and 
(IA18|) is achieved, the energy per particle can be calculated as 

iW-jE.I (l£ + 5 H,p ' (ei ) ■ (A19) 

T \k\<kp T 

The BHF calculation carried out in this paper uses the realistic Argonne V18 (Avl8) [63] nucleon-nucleon interaction 
supplemented with a three-body force of Urbana type, which (for the use in BHF calculations) was reduced to a two- 
body density-dependent force by averaging over the third nucleon in the medium [68| . This three-body force contains 
two parameters that are fixed by requiring that the BHF calculation reproduces the energy and saturation density of 
symmetric nuclear matter. We note that more microscopically based three-body forces without adjustable parameters 
have been recently constructed (see Refs. [69l - l7l| for a recent analysis of the use of three-body forces in nuclear and 
neutron matter). We note also that the Avl8 interaction contains terms that break explicitly isospin symmetry. 
Therefore, in principle, we should consider also odd powers of y in the expansio n fl] ) for the Brueckner calculation. 
However, we have neglected such terms since, as shown by Miither et al. in Ref. |72| . the effects of isospin symmetry 
breaking on the symmetry energy are quite weak (less than 0.5 MeV for a wide range of NN interactions). 
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